Dimensional Study of the Caging Order Parameter at the Glass Transition 
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The glass problem is notoriously hard and controversial. Even at the mean-field level, little is 
agreed about how a fluid turns sluggish while exhibiting but unremarkable structural changes. It 
is clear, however, that the process involves self-caging, which provides an order parameter for the 
transition. It is also broadly assumed that this cage should have a Gaussian shape in the mean-field 
limit. Here we show that this ansatz does not hold. By performing simulations as a function of 
spatial dimension d, we find the cage to keep a non-trivial form. Quantitative mean-field descriptions 
of the glass transition, such as mode-coupling theory, density functional theory, and replica theory, 
all miss this crucial element. Although the mean-field random first-order transition scenario of the 
glass transition is here qualitatively supported and non-mean-field corrections are found to remain 
small on decreasing d, reconsideration of its implementation is needed for it to result in a coherent 
description of experimental observations. 



Introduction - If crystallization is avoided, slowly 
compressed (or supercooled) fluids eventually form a 
glass. They become non-ergodic when their structural 
relaxation timescale r Q (or inverse diffusivity 1/D) gets 
larger than the annealing time. A variety of competing 
descriptions propose to explain this seemingly straight- 
forward process but existing experimental and nu- 
merical results do not allow to unambiguously discrimi- 
nate between them. Yet consensus has recently emerged 
that a growing dynamical length scale is associated with 
the transition, which some have argued results in a crit- 
ical phenomenon of a new kind Q . Based on this devel- 
opment, it seems natural to rephrase the problem start- 
ing from a mean-field (MF) theory, in which correla- 
tions are neglected, and to add correlations progressively 
using renormalization group techniques. Unfortunately, 
even identifying what should be the MF microscopic phe- 
nomenology is fraught with contention. One of the com- 
mon MF framework is the random first-order transition 
(RFOT) theory, which stems from the exact solution of 
a class of abstract spin glass models whose phenomenol- 
ogy is remarkably similar to that of structural glass form- 
ers 0, 5| ■ Inspired by this analogy, reasonable predictions 
have been obtained for realistic models [Bl-Q ■ Yet despite 
these advances, the foundations of the RFOT scenario are 
insufficiently robust for it to be widely accepted as the 
MF theory of glasses, leaving ample room for criticism 
and alternative formulations IlL 

Briefly, the RFOT scenario states that an unavoidable 
crgodicity breaking occurs at finite pressures and tem- 
peratures, whatever the annealing rate. In the MF ap- 
proximation, which is assumed to hold when d —> oo, r a 
diverges at a dynamical transition associated with self- 
caging at which phase space breaks up into pure states. 
In finite dimensions, the growth of r Q at the dynami- 
cal transition is limited, because nuclcation of one glassy 



state from another may be possible up to the Kauzmann 
transition, beyond which only the lowest free-energy state 
prevails [4|. But because nucleation gets strongly sup- 
pressed with d, the dynamical transition is thought to 
dominate the slowdown in high d. 

Like for critical phenomena, RFOT's MF predictions 
are expected to be more accurate above an upper criti- 
cal dimension, shown to be d u = 8 [Tl], 12\- On finite- 
dimensional lattices, a renormalization group analysis 
has even found a RFOT-like fixed point Ijj], although 
in certain cases it disappears in low d [14|. Given that 
these analyses are restricted to abstract models, a key 
question is whether RFOT provides reliable quantitative 
predictions for a realistic particle model, at least in large 
d. This pro gra m was initiated soon after the theory's 



formulation [151 ] . but is no simple task. The key diffi- 



culty is that, while for spin glass models the order pa- 
rameter for the glass transition is the relatively straight- 
forward local Edwards- Anderson overlap [3] , for particle 
systems the caging order parameter is a non-trivial func- 
tion of space, the so-called non-ergodic parameter or its 
Fourier transform, the van Hove function G s (r, t) [lj| . 
The RFOT formulation for particle-based systems re- 
quires a set of integral equations to describe G s (r,t), 
which is only achieved under poorly controlled approxi- 
mations and results in non-equivalent treatments in finite 
d. Most of these formulations, such as density-functional 
theory (DFT) and replica theory (RT) 0, Q, can 
be extended to d — > oo, so if their underlying approxima- 
tions were truly MF in nature they should provide equiv- 
alent and accurate predictions in that limit. Although 
mode-coupling theory (MCT) was developed indepen- 
dently from RFOT [I6| , many have suggested that a dy- 
namical description of the RFOT scenario should result 
in MCT- like equations I5|, ll7( . This observation raised 
the question whether MCT should converge to the cor- 
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FIG. 1. Power-law fits (lines) to the vanishing diffusivity 
improve with d=4-9 (different symbols), spanning over three 
decades of D. (top inset) The resulting numerical <pa values 
(including the d = 10—12 results from Ref. [2^1 ) ar 6j however, 
significantly different from the MCT results for the dynamical 
transition (short-dash line), (bottom inset) The results for 7 
also disagree with MCT predictions. 



rect mean- field description in high d [15[. MCT's results 
for the glass transition were, however, recently found to 
be not only asymptotically divergent from DFT's and 
RTs [3 El, but t0 S et 

increasingly unphysical even in 
relatively low dimensions (20| . In order to evaluate the 
MF scenario for the glass transition, we here stringently 
test these theories against simulation results as a func- 
tion of d, emphasizing the evolution of the caging order 
parameter. 

Results - Hard sphere fluids are the simplest glass 
formers with which to compare theoretical prediction, 
because their structure gets increasingly trivial with 
d 0, H2]. As a first control, we test MCT's power-law 
scaling form for the vanishing diffusivity D ~ (tp — </Sd) 7 
when the fluid packing fraction ip approaches the dynami- 
cal transition at tp^. This form fits high-d results well, ex- 
cept for the more sluggish systems in d = 4 (Fig. [T]) [23j . 
and the ipd value s ag ree with those obtained from a dif- 
ferent procedure [20(. MCT predicts values for tp^ and 
7, however, that are inconsistent with the numerical re- 
sults, while RT predicts a Kauzmann transition bound 
> consistent with the numerical data l2(il l. but 



results for 7 arc still missing [24 1. 

The caging order parameter at the dynamical transi- 
tion should offer a clearer picture of what is happening. 
In the high-d MF limit, the van Hove function is argued 
to be Gaussian based on the multiplicity of "uncorre- 
lated" caging neighbors and the central-limit theorem. 
In practice, most implementations of both RT 0, 0] and 
DFT UK simply assume a Gaussian form; while MCT 
predictions have been understood as faulty partly be- 
cause they do not tend toward one [lij]- Belief in the 




FIG. 2. (a) The MSD in d = 6 for increasing </>=0.1453- 
0.1720 illustrates the developing caging regime (dashed line), 
intermediate between the ballistic (thick line) and the diffu- 
sive (solid line) regimes, (b) The MSD for isodiffusive states in 
d=4-8 identifies the caging mid-point time tmp (large dot). 
The power law fitting parameters for the caging regime in 
d = 6 from (a) are used in (c) to extract the plateau height 
at the dynamical transition, when f = 0. (d) The plateau 
height (solid line) is consistent with RT's lower bound at ip-K 
(long-dash line) and significantly different from the MCT pre- 
dictions (short-dash line). 



Gaussian form is so anchored that sustained deviations 
from it were quickly interpreted as dynamical heterogene- 
ity absent from the MF picture 

USUI- Directly mea- 
suring the van Hove function for hard spheres at ipd, i-C, 
once diffusion is fully suppressed, is, however, impossi- 
ble. In low d, the dynamical transition from the RFOT 
scenario is avoided, which transforms the arrest into a 
dynamical crossover and blurs its properties. Although 
increasing dimensionality resolves this ambiguity, reach- 
ing equilibrated configurations near ipa is challenging. 
Annealing more slowly than r Q limits the numerically 
accessible r Q to those within a few orders of magnitude 
from the collision time. 

To circumvent this difficulty, we examine the system- 
atic development of the caging regime of the mean square 
displacement (MSD) when approaching the dynamical 
arrest as d increases. Near ipa the MSD develops an in- 
flection between the ballistic and the diffusive regimes, 
which should plateau at full caging (Fig. [2jt). On a 
purely phenomcnological basis, we describe this interme- 
diate regime by a power law 



(Ar 2 (t)) = 2da(tp)t c ^ , 



(1) 



whose subdiffusive exponent < 1 decreases with 
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(a) 
d=6 



H>—Q. 1453-0. 1720 




FIG. 3. (a) The evolution of the van Hove function with 
packing fraction at tmp in d = 6 shows that the fat expo- 
nential tail (short-dash line) steadily grows at the expense of 
the Gaussian regime (long-dash line), (b) The isoconfigura- 
tional results for four randomly chosen particles (symbols) at 
tp = 0.1720 in d — 6 indicate that the individual cages as well 
as the average cage (thick line) are non-Gaussian, (c) The 
isodiffusivity indicates that the fat tail remains undiminished 
for all the dimensions studied, and grows increasingly differ- 
ent from the MCT results (long-dash lines), given here for 
d=4, 6, and 8. 



increasing tp. Under the reasonable assumption that 
C(v) ~ ^ for tp — > (pa, extrapolating the paramet- 
ric plot a(C) to the limit ( —> gives the cage size 
A = a(( — > 0) (Fig. f2b). This measure is found to 
remain essentially constant 2dA m 0.027(4) over the d 
range considered. The MCT description of the interme- 
diate regime is in qualitative, although not quantitative, 
agreement, and RT's lower bound from caging at tpn is 
respected (Fig. [5H). 

Considering the caging order parameter, rather than 
the cage size, more directly probes glass formation. We 
consider the evolution of G s (r, tmp) at the logarithmic 
mid-point tmp of the intermediate caging regime. This 
choice has the advantage of correctly extrapolating to the 
full caging limit at the dynamical arrest, while staying 




FIG. 4. The isodiffusive for d =4-8 (main) and d = 6 with 
increasing tp (inset) non-Gaussian parameter a.2(t). The peak 
at r a decays strongly with d, an indication that dynamical 
heterogeneity is increasingly suppressed, while the caging be- 
havior near tmp is robust. 



well clear of r Q (Fig. [5Jd). The results surprisingly indi- 
cate that the cage shape does not tend toward a Gaus- 
sian (Fig. [3]). The Gaussian regime in fact shrinks to 
smaller r with increasing d and tp, leaving instead a re- 
markably fat tail. This result contrasts with the RT and 
DFT assumptions, and markedly differs from the MCT 
predictions, whose discrepancy grows worse with d. 

One may wonder if this pronounced deviation from the 
expected Gaussian behavior is due to dynamical hetero- 
geneity, and to the growth of an associated dynamical 
length scale, as has been found in low d @, O [27l429j ] . 
The RFOT scenario and a description based on dynami- 
30l . 31 1 both suggest that such a length 



cal facilitation [1( 
scale should be present in all d, but like other critical 
lengths, the distance (tp — </?d)/Vd from the critical point 
over which it is felt is expected to shrink with d [ll , [l2| . 
The impact of dynamical heterogeneity should thus effec- 
tively disappear with increasing d. The van Hove func- 
tion unambiguously resolves the two processes. We first 
perform an isoconfigurational study, in which a same ini- 
tial configuration is randomly assigned series of different 
random initial momenta [27j , indicating that the individ- 
ual particle cages are not Gaussian either (Fig. Et") - If the 
non-Gaussianity arose from an heterogeneity of the local 
relaxation on the tmp scale, then one would expect the 
individual cages to be Gaussian, which they are not. We 
then consider the non-Gaussian parameter a 2 (t), which 
is the kurtosis of G s (r,t) (Fig. g}. Although a 2 (t) de- 
creases with d for all time regimes, for isodiffusive sys- 
tems the change is much more pronounced at the peak 
of CK2 (t) near T a , where dynamical heterogeneity is max- 
imal, than at the caging mid-point tmp- These results 
therefore support the numerical evidence that G s (r, tmp) 
remains no n- Gaussian in the MF large d caging regime. 



4 



Discussion - Our work clarifies the MF scenario of 
the glass transition and establishes mileposts for assess- 
ing current and future theoretical descriptions of the phe- 
nomenon. The results suggest that the RFOT scenario 
qualitatively describes high-dimensional hard spheres, 
and that non-MF corrections remain small upon de- 
creasing dimension, even below d u . When d increases, 
the power- law divergence of 1 / D near tpd is clearly visi- 
ble; the associated dynamical heterogeneities around r a 
meanwhile decrease, making the glass transition a local 
caging problem describable by MF theory. Yet, contrary 
to common belief, local caging does not lead to a simple 
Gaussian caging order parameter. 

We hnd a smooth dimensional dependence of the struc- 
tural and dynamical properties, which is consistent with 
what is found in the dynamical facilitation scenario [32| . 
One might thus wonder whether a smooth d dependence 
disagrees with the RFOT picture, which is based on an 
underlying critical phenomenon with an associated upper 
critical dimension d u = 8 11, 12f . Yet as in standard crit- 



TABLE I. Numerical properties extracted from simulations. 



d <p d 



A 



ical phenomena, the Landau-Ginzburg criterion indicates 
that one should be extremely close to the critical point to 
see deviations from the MF predictions, even below d u . 
Recent quantitative computations show that the regime 
where non-MF corrections are present is quite hard to 
access using numerical simulations in d = 3 0, [24[ , and 
it is reasonable to argue that it should be even harder to 
reach with increasing d. The existence of d u is therefore 
expected to be undetectable unless one is much closer to 
the critical point than we are here, and our results thus 
remain qualitatively consistent with the RFOT scenario. 

Despite this qualitative agreement, we show that all 
the concrete implementations of RFOT theory struggle to 
describe the high d regime, although it should be the eas- 
iest. We find that a broad scope of MCT predictions are 
defective: the predictions for i^d, the exponent 7, and the 
cage shape are not only wrong, but worsen with increas- 
ing d. These results reveal the inadequacy of standard 
MCT as a MF description, challenge some of the deep- 
seated assumptions about glass formation, and strongly 
call for a revised formulation of a dynamical theory of 
the RFOT [H. At the same time, DFT and RT assume 
from the very beginning a Gaussian form for the cage, 
which is invalidated by our results. Even if some of the 
RT results seem consistent with our computations, the 
theory should also be revised to understand the extent 
to which a non-Gaussian caging order parameter affects 
its predictions. 

If the RFOT formulation is indeed correct, we expect 
that theoretical reconsiderations will lead to a resolu- 
tion of the discrepancy between the MCT and RT/DFT 
predictions [3, [l9j], and ultimately to a consistent de- 
scription of both the statics and dynamics of glass for- 
mers. Given that non-MF corrections seem small even 
for d < d u , it is even possible that such a complete MF 
theory could perform quite well in experimentally rele- 



4 0.4065 2.38(6) 0.027(2) 

5 0.2700 2.34(4) 0.027(2) 

6 0.1732 2.25(4) 0.030(4) 

7 0.1081 2.22(3) 0.027(2) 

8 0.06583 2.12(3) 0.029(5) 

9 0.03938 2.10(3) 0.026(5) 



vant dimensions. It would therefore provide a productive 
starting point for a more refined renormalization group 
analysis [13| that takes into account the role of fluctua- 
tions below d u - It would also be interesting to evaluate 
the high-dimensional robustness of descriptions based on 
dynamical facilitation [10|, |30|, |3l| . This program seems 
a promising route for obtaining a more robust and less 
controversial theory of glasses. 



MATERIALS AND METHODS 
Numerical Simulations 

Event-driven molecular dynamics simulations of 8000 
hard spheres in dimensions 4 < d < 9 are performed 
under periodic boundary conditions (2(| 33 1. Because 
crystallization in high d is strongly suppressed, ac- 
cess to deeply supersaturated starting configurations can 
be achieved via the slow compression of a low-density 
fluid [1^, Between 4 and 8 independent configura- 
tions are equilibrated for each packing fraction. Simula- 
tions are then run at constant unit temperature kgT for 
times at least as large as 10/(2dD), where time t is ex- 
pressed in units of (ma 2 /ksT) for particles of unit mass 
m and unit diameter a. Even in d = 9 for the densest 
system studied, the box side is kept greater than 2a. Be- 
cause the static and dynamical correlation lengths shrink 
with d, these system sizes avoid significant finite-size ef- 
fects, as discussed in 12011 . 



The average mean-square displacement (MSD) 



A' 



(Ar^HTV-^dr^-r^O)! 2 



(2) 



is obtained from equilibrated starting configurations. At 
times shorter than the collision time, MSD displays a 
ballistic regime (Ar 2 (t)) = dt 2 , and at long times it has 
a diffusive regime (Ar 2 (t)) = 2dDt. 

Fitting these numerically determined D to the power- 
law form D ~ (<p— (fid) 7 is reasonably good for D < 0.005, 
and improves with increasing d. The resulting values 
of 7 and (fid are reported in table 1 and Fig. [TJ For 
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d > 5 [23[, the full accessible dynamical range stud- 
ied is used, spanning up to three D decades. Isodiffu- 
sive comparisons are made for systems whose <p gives 
2dD = 0.0008(3). Because 7 differs relatively little 
over the d range studied, choosing isodiffusive systems 
is roughly equivalent to keeping the distance to the dy- 
namical transition (tp — cp^)/^ constant. Note that small 
differences in 2dD can affect some measures, such as the 
non-smooth evolution with d of the peak near r Q in Fig.|U 
The cage is described by the self part of the van Hove 
function 

1 N 

G s (r;t) = -J2{S(\r t (t)-r, l (0)\-r)) , (3) 

i=l 

which in the ballistic and diffusive regimes is well- 
approximated by a pure Gaussian. The logarithmic 
caging mid-point tmp is chosen at the mid-time on a 
logarithmic scale, intermediate between the ballistic and 
the diffusive extrapolations of the MSD. The isoconfig- 
urational study of the cage was repeated for 1000 dif- 
ferent random initial velocity distributions, in order to 
obtain a good statistics on the individual cages. The 
non-Gaussian character of this distribution is canonically 
described by its kurtosis, or non-Gaussian parameter, 



a 2 (t) 



<Ar 4 (*)> 



d + 2 (Ar 2 (t)y 



1 . 



(4) 



Mode-Coupling Theory 



Replica Theory 

The best replica scheme for studying hard spheres is 
the small cage expansion (tJ. It consists of taking the 
lowest-order expansion of the replica theory free energy in 
the cage size A, as given in Ref. @, Section VII, Eq.(73)]. 
This approximation gives reliable results at high density 
near jamming, where A is small. Unfortunately, using 
this scheme the dynamical transition, which corresponds 
to the point where the self-consistent solution for states 
with cage A vanishes, is not found. Because of the crudc- 
ness of this approximation, the equation for A has indeed 
a solution A*(ip) at all (p 0, Eq. (74)] . Taking into account 
higher orders in the small A expansion is only possible 
in the limit d — > 00, leading to an asymptotic prediction 
for the dynamical transition 0. However, this asymp- 
totic limit is reached only for extremely high dimensions, 
d > 50, that are not accessible to simulations [2(|. Other 
replica schemes are available 0, but they do not give 
good quantitative results. In summary, for the moment 
replica theory does not give reliable predictions in the 
regime that is relevant for the present study, namely 
ip ~ ipd and low d. Using the lowest order expansion in 
A Section VII, Eqs.(73)], we can nonetheless obtain 
the cage radius at the Kauzmann transition tpx > <p&, 
which provides a lower bound A*(cpd) > A*((pn) for the 
cage size. 

We thank J. Kurchan and R. Schilling for stimulating 
discussions. P.C. acknowledges NSF support No. DMR- 
1055586. 



The MCT analysis follows the approach of Refs. [18, 
[l9l | , using the Percus-Yevick (PY) structure factor calcu- 
lated iteratively with a numerical Hankel transformation 
of order d/2 — 1. The agreement between the structural 
PY prediction and the numerical results improves with 
d, so the MCT predictions are not expected to depend 
on this choice of input. Using the hypcrnctted chain 
(HNC) input for the structure factor only worsens the 
agreement with simulations. The long time limit of the 
self part of the van Hove function G s (r;t —> 00) is cal- 
culated by Fourier transforming the MCT solution of 
the self part of the non-ergodic parameter f s (q). The 
plateau height at the long time limit of the MSD at 
caging 2dA = lim f ^ 00 (Ar 2 (t)) is calculated through a 
small wave number analysis of the MCT equation, 



dq 



q 2 c(q) 2 S(q)f(q)f s (q), 



(5) 



where f(q) is the collective part of the non-ergodic pa- 
rameter, for a given static structure factor S(q) and direct 
correlation function c{q). In order to check the consis- 
tency of the numerical calculation, we also obtained the 
plateau height through 2dA = J df r 2 G s (r; t —> 00) and 
found the relative error to be smaller than 1%. 
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